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Abstract. This paper derives a general procedure for the numerical solution of 
the Lindblad equations that govern the coherences arising from multicoloured light 
interacting with a multilevel system. A systematic approach to finding the conservative 
and dissipative terms is derived and applied to the laser cooling of gallium. An 
improved numerical method is developed to solve the time-dependent master equation 
and results are presented for transient cooling processes. The method is significantly 
more robust, efficient and accurate than the standard method and can be applied to 
a broad range of atomic and molecular systems. Radiation pressure forces and the 
formation of dynamic dark-states are studied in the gallium isotope ^^Ga. 



1. Introduction 

Many recipes [1-5] now exist in the literature for ultracold molecules starting from 
a small selection of basic ingredients: trapped, ultracold atoms [6]. The variety of 
molecular species is limited by the fact that, to date, only a small number of elements 
have the requisite electronic structure for direct laser cooling, although sympathetic 
cooling of molecules with ultracold buffer gases may prove an effective alternative 
method [7]. In order to increase the repertoire of ultracold molecules available, 
techniques must be developed to cool the large swathes of the periodic table currently 
inaccessible to laser coohng. 

The p-block elements form a particularly important part of the periodic table, but 
the only such atoms that have been cooled and trapped thus far are the noble gases, 
paradoxically the least reactive of all elements. There has been, however, some evidence 
of laser cooling and manipulation of elements within Group 13: a beam of aluminium 
atoms [8,9], has been observed to narrow under laser irradiation and a deceleration force 
has been demonstrated in gallium [10] and indium atoms [11] with recent evidence of 
sub-Doppler transverse cooling in an indium atomic beam [12]. 

In general two problems need to be addressed for cooling in the p-block. The 
first is the tendency of all these elements to absorb light in the UV rather than the 
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visible region when in their ground electronic states. Stable continuous-wave lasers 
that emit in the UV tend to be complex to operate and generate relatively low power. 
Sometimes it is possible to find a more convenient, laser- accessible, transition starting 
from a metastable state, which has been the cooling strategy in the noble gases [13,14] 
and in the successful Group 13 studies conducted so far. The second is the p-shell 
itself, with the associated fine structure. The energy splitting due to spin-orbit effects, 
is many times greater than the hyperfine splitting and necessitates more than one laser 
(frequency). The usual strategy of repumping can lead to coherent dark states that 
compress the cooling process [15]. 

Group 13 elements fortunately are untroubled by the first problem, which is more 
an issue with laser technology rather than fundamental physics, having absorption 
frequencies in the visible and near-UV (with the exception of boron). The (ns^,np) 
configuration has a simple doublet ground state, thus the presence of spin-orbit coupling 
in the affords the best opportunity to address the practical issue of cooling in the 
presence of complex electronic structure. In our first paper [16], we conducted time- 
dependent calculations on the interaction of ^^Ga atoms with counter-propagating laser 
beams, but the method involved was computationally inefficient, required the derivation 
of the equations of motion manually and did not easily allow the exploration of optimal 
experimental parameters, the information experimentalists would find most useful. In 
this paper, a more sophisticated and elegant method has been developed to quickly and 
efficiently evaluate possible cooling schemes within these complex structures. 

2. Mathematical Model 

We employ a standard semiclassical treatment throughout, and in this section we briefly 
introduce the notations, definitions and approximations that we employ in this work. 
The notation is important in forming, in a systematic manner, the structure of the 
differential equations both mathematically and computationally. 

The state of the atom is expressed in a finite set of Na eigenstates of the hyperfine 
Hamiltonian (H^) which we denote by the single index, p G {1,2, . . . , N^}. This 
abbreviation stands for the ensemble of the hyperfine sub-level quantum numbers. 
Expanding this notation, we take, ctp to denote the collective label for the fine- 
structure level and uncoupled nuclear spin. That is, ap — {up, Lp, Sp, Jp, Ip} where 
Up is the electron principal quantum number, {Lp, Sp) are the orbital and spin angular 
momenta with the resultant coupling Jp, and Ip is the nuclear spin. The basis functions 
\p) — \ap, Fp, Mp) are the orthonormalized eigenstates with eigenvalues: Ep (in general 
not unique). 

H^\p)-ep\p) (1) 

For convenience we define the corresponding angular frequencies: ficup = £p. 

The atoms are moving at non-relativistic speeds at all times with the centre-of-mass 
of the atom moving with a velocity v in the laboratory frame, and thus the Galilean 
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transformation can be employed. The internal state of the atom is expressed by the 
(variation of constants) expansion: 

Na 

I^W) = E<Wlp) (2) 

Then the density operator, p= \^{t)){^{t)\, has the matrix representation in the time- 
stationary basis set: 

Pi,{t)^{p\p\m = c^{t)cy{t) (3) 

The multicoloured lasers are represented by a superposition of classical 
monochromatic fields: 

Nl 

-^(^) ^) = ^3 cos(fej • r - Ujt + Lpj) (4) 
j=i 

with the polarizations, wavevectors, freqencies and phases denoted by the usual symbols. 
The superpositions implicitly allow for any state of polarization. Then the time- varying 
perturbation can be written: 

H^''{t) = -d-E{r,t) (5) 

where d is the (electronic) dipolc-moment operator. 

The evolution (Liouville) equation, in the absence of dissipation, takes the form: 

|p^(t) = -{t/h) [h^ + H^\t),p^{t)] (6) 

The hyperfine Hamiltonian is diagonal in this basis, and in the absence of the 
couphng, H"^^ — 0, we have the solution. 

p^(t) = e-'"^'p^(0)e'"^' (7) 

In other words, the density matrix in the interaction representation is stationary 
with the form: p^/(0) = (0)Cp, *(0). The dipole approximation is vahd for these 
frequencies, and thus for an electron with coordinate rg with respect to the centre-of- 
mass of its atom, and thus with a position r -|- re in the laboratory frame, we have: 
kj ■ {r + Ve) ^ kj • r — kj -vt + cl). Then the couphngs can be written in the form: 

H^it) = J2 ^^j,PP' cos{kj ■ vt - Ujt + ip'j) (8) 

which reveals the Doppler shift in a simple manner. The Rabi frequencies determine 
the selection rules and are defined according to: 

^j,PP' = -^Ej ■ {p\d\p') (9) 



It is convenient (analytically) to transform to the interaction representation, 



Pi,, ^ e^^'^p^t)^^'^ (10) 



and with V{t) = e^^^^H^^{t)e we have the evolution equation: 

y{t) = -{i/n)[v{t),p\t)\ (11) 
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At this point it is conventional, though not necessary, to make the rotating-wave 
approximation and discard non-resonant transitions. That is the exponents of the 
commutator involve the Doppler-shifted detunings 

^tpp' =^P- ^P' ±kj-v^ Uj (12) 
where ± corresponds to absorption (+) or emission (— ). Once this is done, a 
further unitary transformation to remove the time-dependence (oscillations) from the 
interaction is performed, such that: H = U(t)V(t)U''{t). This is only possible when the 
rotating-wave approximation is enforced, and means that the time-integration is now 
relatively straightforward. 

The inclusion of the spontaneous emission terms can be modelled by the Lindblad 
correction which is written in the form: 

^^pit) = ~iz/h)[H,pit)]+T.p (13) 

where the dyadic (tensor) for spontaneous emission is denoted by F and expresses the 
decoherence dynamics of the system. 



3. Outline of computational method 

The master equation (13) provides the dynamics of the model, within the rotating-wave 
approximation. This equation is used to study the transient and steady-state behaviour 
of the system. However, naive numerical approaches, such as the Euler method, can 
lead to unstable integration. We propose the following combination of factorisation and 
modificd-Euler method as the integration algorithm. 

• Interaction matrix elements 

The set of atomic hyperfine states { \p)} is defined in terms of the quantum numbers, 
and energies. The selection rules for absorption and stimulated emission are 
enforced by the computer program based on this data. Then, for each (allowed) 
transition, the multicolour detunings are specified and thereby the non-resonant 
couplings (counter-rotating terms) identified and removed. 

• Coupling terms 

The reduced matrix elements are assumed as known, so that one simple uses the 
Wigner-Eckart theorem [17] to evaluate the coupling coefficients in terms of the 
reduced matrix element. 

For convenience we introduce the (reduced) Rabi frequency: 

^j,PP' ^ ^Md\\p')\ (14) 

• Dissipation 

The (partial) lifetimes of the excited states are assumed to be known. That is the 
spontaneous decay process from any upper state to any lower state: \p) \p'), is 
given as: 

2^-' = 3 (2F, + l)^c3 
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So for each of the dipole allowed transition the relative strength of stimulated and 
spontaneous emission is given by the saturation parameters for each transition: 

Gj,pp' = '^^hp'Hp' (16) 
Numerical integration 

Consider (13). Given the initial state, p(0), we step forward in time to t — h, using 
the intermediate helf-step as follows: 



.1 



hH 



P(0) + IhT . p(0) 



The method is still accurate (only) to first-order in h, nonetheless it is much superior 
in stability to the standard Euler method. There is the computational cost of the 
exponential operators (though these have to be calculated only once within the 
RWA), but in terms of numerical accuracy and efficiency the method is a significant 
improvement, as will be shown in section 4. The steps are repeated in the usual 
manner to produce the entire evolution over time. The Rabi oscillations are rapid 
at the beginning followed by slow damping so that small time steps ( Q/i <^ 1 and 
7/1 <^ 1) are required to ensure that the integration proceeds correctly. However, at 
longer times, (7^ » 1) as the steady state approaches, the characteristic oscillations 
are slower, and larger time steps can be taken. 

Spontaneous emission 

As the spontaneous emission tensor, F, depends on the density matrix, p, it must 
be calculated at each time step; this involves solving the equations for the atomic 
density matrix elements pki = PaaFaMa,attiMt = {aaFaMa\p\abFbMb) . To solve the 
equations each possible pair of ground and excited-state sublevels are considered in 
turn, and certain selection rules are enforced: 

AF = 0, ±1 (but not F = ^ F' = 0) 

AM = 0,±1 (19) 
Mg- Mg= even 

Mg-M'g^M,- Ml (20) 

Each element of the matrix, F, can be described by one of the following four 
equations [18]: 



= -(7aeiFei + la,^F,^){ae^Fe,MeM(^e2Fe2Me2) (21a) 

{aeFeMe\V p\agFgMg) = --f^^FMeFeMe\p\agFgMg) (216) 
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Oej ,Oe2 ,Fei ,Fe^ ,Mei 

X {ae,Fe,Me,\p\(^e,Fe,Me,) (21c) 

{a,FgM'g\Tp\agFgM,) 

= 5] {FgMgM'g\A\FeMeM'^){aeFeM'MaeFeMe) {21d) 

where 

X ^ {Fg,MgMFeM{Fg,Mg,lq\Fe,M,,) (22a) 

g=0,±l 

(F^M.M^I AjFeMeM^) 

^'i^a.F.^a.F, E (i^.M^l^lFeM^) (F,M,lg|FeMe) (226) 
?=0,±1 

(F,M,lg|FeMe) = (-l)^-^+^V2Fe + 1 ( J (22c) 

The Wigner 3j symbols (22c) are calculated using the MATLAB code wignerSj.m 
[19]. 

Using this method, the time-evolution and steady-state of the density matrix can 
be studied. As well as containing information about the populations of each state, the 
density matrix is also used to calculate the dipole radiation force on the atoms. This 
force, F, is due to the interaction of the laser field with the induced atomic dipole 
moment and is given (in one dimension) by; 

F{z,v,t) = -^U{z,v,t) (23) 

where U{z,v,t) = —tY{pd)E{z,t), E[z,t) is the electric field of the lasers along the 
propagation direction and d is the atomic electric dipole operator (in the atom frame). 

4. Analysis of Method 

To analyse the accuracy and stability of our method, it was compared to the basic 
Euler method [20] . As the exact solution of the density matrix equations for a two- level 
atom, with ground state, b, and excited state, a, is well known [21], this system was 
used to compare the two methods. Consider an atom at rest (f = 0), initially with all 
population in the ground state, Paa — 0. The population of the excited state, Paa, can 
be calculated as 

" (2G + 8) ~ cosAt+ ^sinAt exp(-|7i)| (24) 

where G is the saturation parameter, G — 20^^/7^, with 7 the spontaneous decay rate 
and A = JV2G- 1. 
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Time, yt 

Figure 1. Population paa for tlic excited state of two-level atom, in which Paa{0) — 0, 
G — 10, S = 0, J = 1, h = 0.1, at velocity fc-y/7 = 0. Figure shows the results of the 
numerical integration method presented in this work (dashed) , compared to the Euler 
method (dot-dashed) and the exact solution (solid). 

Figure 1 shows the time-evolution of the excited-state population paa, calculated 
on resonance {kv = 5 = 0) by both the Euler method and the method presented in 
this work, compared with the exact solution given by (24). It is clear that even for 
this simple system the new method converges significantly more rapidly than the Euler 
method. This slow convergence of the Euler method becomes even more marked as 7 
is reduced. At 7 = (pure Rabi oscillation), and in calculations far off-resonance, the 
Euler method fails to converge unless step sizes are extremely small {h <0.001) whereas 
the method presented is highly stable, even at step sizes over 100 times larger. 

Comparisons for more complex systems were not possible as the extremely small 
step-sizes needed to counteract the instability of the Euler method meant calculations 
could not be carried out in realistic timescales. In contrast, the method presented here 
remained stable over a wide range of detunings and step-sizes. It is clear therefore 
that the new method provides a marked improvement on accuracy, stability and 
computational cost. 

Previous calculations [16] on a similar system to that outlined in this paper 
depended on the construction, by hand, of a set of n x n density matrix equations, 
where n is the number of individual magnetic sublevels in the cooling scheme. For small 
systems this is a relatively simple task but as the number of magnetic sublevels increases, 
so does the number of equations needed. Soon it becomes hugely time-consuming to 
construct and solve these equations, with the potential for human error also increasing. 
The most difficult part of the master equation (13) to calculate is the section containing 
the spontaneous emission terms, r(p), which often has to take into account many 
different decay channels, alongside their associated Clebsch-Gordan coefficients. The 
benefit of the method presented here is that, with the only input being a small number 
of quantum numbers, these terms are calculated automatically within the numerical 
code, and placed in a matrix. This abolishes the need for time-consuming construction 
of equations or matrices by hand. 
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5. Application to a A-system 

To show the possibihties of this new method, in this paper the method presented is 
apphed to a simple multi-level system. Consider now a A-system involving the lower 
manifolds |g), |G), and the excited state |e), as shown in figure 2. This system is 
assumed to correspond to the — )■ transition in a Group 13 atom with nuclear 
spin 1=0, specifically ^^Ga, a beta emitter and putative radiopharmaceutical for PET 
imaging which can be produced with a biomedical cyclotron [22]. This scheme could 
also be applied to a transition to the first hyperfine excited level of a Group 1 atom 
with nuclear spin 1=1 (eg ®Li). 



For the example presented, two laser frequencies are considered, ooi and 002, red- 
detuned by 61 and 62 respectively. One laser beam corresponds to the Fg = 1/2 — >■ 
Fe = 1/2 transition, with the other used to pump Fq = 3/2— )'Fe = l/2. The scheme 
considered here assumes the use of counter-propagating circularly-polarized laser fields, 
with a cr^-polarized beam driving the Mp — )■ Mf + 1 transition, and a a~ laser driving 
Mp ^ Mp- 1. 

For simplicity it is assumed that the radiative widths for the two transitions are 
approximately equal, 71 ~ 72 = 7, that ki ^ k2 = k, and that the saturation parameters, 
Gi = G2 = 1- Tests with a range of parameters show that these assumptions have little 
effect on the results of the calculations. 

5.1. Population Distribution 

The population distribution of the system is investigated over a range of atomic 
velocities, along with the dipole radiation force, for both the time- dependent and steady 
state cases. The steady state result is obtained numerically, for any initial conditions, 
by time-dependently solving the master equation at very long timescales, until the 
populations converge. Figure 3(a) shows the steady-state ground-state population 
distribution for a 2-(4)-2 multilevel atom, as shown in figure 2, with laser detunings 
61 = —37 and 62 = —27. 



e-1/2 ei/2 




g-1/2 gi/2 



Figure 2. Schematic energy level diagram for 2-(4)-2 multilevel atom, corresponding 
to the ^Pi/2, ^^3/2 ^ transition in Group 13 atoms with nuclear spin 1=0. 

corresponds to the Rabi frequency, uj is the laser frequency and S is the laser detuning. 




Figure 3. Ground-state populations g±i/2 (solid), G-i-1/2 (dashed) and G-i-3/2 (dot- 
dashed) of a 2-(4)-2 multilevel atom (figure 2) in the steady state as a function of 
velocity v — v^. Radiative widths are equal, 71 = 72 = 1, as are laser intensities, 
Gi = G2 = 1, while detunings are Si = —87, 62 = —'2"f. The figure shows sharp 
two-photon resonances at velocities fcu/7 = ±0.5 and kv/^ — 0. Figure (b) shows 
an enlarged view of (a), with G = 1, but as the population distribution is highly 
symmetrical the states g-1/2, G_i/2 and G_3/2 are omitted from (b) for clarity, (c) 
also shows an enlarged view, in this case with G — 10, and at this higher laser intensity 
broadening of the two-photon population peaks is observed. 



In figure 3(a), it can be seen that as the velocity of the atom nears resonance 
velocities {kv/'y = ±2, kv/'-f = ±3) there is a transfer in populations giving rise to a 
broad peak, as expected. The system possesses a left-right symmetry, with a change 
in the sign of the velocity being equivalent to a change in the sign of the quantum 
number Mp. Sharp population structures are also observed closer to zero velocity. These 
resonances are due to two-photon processes producing ground-state coherences which 
vary rapidly at low velocities, leading to sharp variations in ground-state populations. 
The positions of these resonances can be predicted simply from the energy conservation 
law [23]. Two-photon transitions within the ground state |G) do not change the energy 
of the atom, so {u2 ± kv) — =F kv) = and the resonance occurs around zero velocity. 
However, in the case of two-photon transitions between states |g) and |G) the energy of 
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the atom is changed; 

{uji ± kv) - {lo2 =F kv) = {loi + 6i) - {002 + 62) (25) 

where Ui is the laser frequency used and 6i is the detuning (as shown in Figure 2). These 
resonances, therefore, occur at velocities kv = ±{Si — 62) /2, in this case at kv = ±0.57. 
Figure 3(b) gives an expanded view of the region in which these population shifts occurs, 
and as figure 3(c) shows, increasing the laser intensity leads to a broadening of these 
peaks. Conversely, a low saturation parameter gives rise to much narrower resonance 
structures. As the population distribution is highly symmetric, for clarity only the levels 
gi/2 5 Gi/2 and G3/2 are shown in the expanded region. 

Steady-state excited populations are shown in figure 4; due to symmetry, the 
population distribution for each excited-state sublevel is identical. Ground-state 
coherences caused by two-photon processes result in the populations falling sharply to 
zero at these resonance points. The presence of these structures at all laser intensities 
clearly demonstrates the destructive atomic interference occurring in the system, with 
the formation of dark states. 



6 

CO 

c 

O 

Q. 
O 

Q- 2 



"1 -5 5 1 

Normalised Velocity, kv / y 

Figure 4. Steady-state excited populations as a function of velocity v = Vz- The 
populations for each of the excited states is identical, with sharp drops to zero 
observed as a result of two-photon resonance processes. These demonstrate destructive 
interference occurring with the formation of dark states; saturation parameter G — 1, 
detunings Si = —87, S2 = —27. 

We now consider the formation of these resonances and dark states by studying 
the time-dependent evolution of the populations. Figure 5 shows how the ground-state 
population distribution varies with time at different velocities. In each of the figures the 
initial conditions are taken as the Boltzmann distribution in a gas of gallium atoms at 
1500K, the typical temperature of a gallium effusive cell. Initially, 34.4% of the atoms 
are in each of the ^Pi/2 states, with 7.8% in each of the ^P3/2 sublevels. At kv/'j = 0, (a), 
a symmetry is observed as expected, resulting in three sets of overlapping populations. 
Due to the two-photon processes and the formation of dark states within the state |G), 
the population in the |g) state decays to zero. As the velocity decreases slightly to 
kv/j = —0.25 (b), the six different populations can be seen to reach their steady-state 
distribution quickly. As the velocity is now non-zero, the symmetry of the system is 
broken, and the behaviour of individual populations can be observed. Figure 5(c) shows 
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Figure 5. Evolution of ground-state populations g_i/2, gi/2 (solid), G_i/2, G1/2 
(dashed) and G_3/2, G3/2 (dot-dashed) at (a) kv/'y ~ 0, (b) kv/'f = —0.25, (c) 
kv/j — —0.5 (two-photon resonance velocity) and (d) kv/^ — —3 (at a resonance 
velocity). At kv/'y — 0, (a), the symmetry of the system dictates that sublevels with 
the same modulus Ai^pj have identical populations. This symmetry is broken for 
non-zero velocities, as shown in (b). At two-photon resonance velocities, as illustrated 
by (c) kv/y = —0.5, markedly different behaviour is observed, particularly in the 
5-1/2 population, which rises to a maximum but falls in every other case. Population 
oscillations also occur. In subfigure (d), taken at a resonance velocity kv — 61, rapid 
changes in population are observed at short timescales (7t < 20). The state g_i/2 
is depopulated quickly, with a resulting gain in population for the G_i/2 and G3/2 
sublevels. Similarly, population in the G_3/2 state decreases, resulting in a gain for the 
gi/2 sublevel. After this time the populations approach their steady-state solution, with 
states with Mp > most highly populated: saturation parameter G — 1, detunings 
Si — —87, S2 — —27, initial ground-state populations; pgg = 0.344, pqq = 0.078 



the population evolution at = —0.5, one of two-photon resonances between the 

two fine structure states. In this figure oscillations are observed, before the populations 
of states G3/2 and g_i/2 rise to a maximum. All other populations tend to zero. This 
behaviour is notably different from that at other velocities. In subfigure (d), the internal 
processes can be seen clearly. At first (7t < 20) population is transferred rapidly from 
the g-1/2 state, resulting in a gain in population, via the excited state e_i/2, for the 
G3/2 and G_i/2 sublevels. Similarly, population in the G_3/2 state decreases, resulting 
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in a gain for the gi/2 sublevel. After this time, the populations approach their steady- 
state distribution, with the gi/2, G1/2 and G3/2 states populated while the other three 
ground-state populations tend to zero. 



5.2. Dipole Radiation Force 

The radiation pressure force on an atom, using (23), can be easily calculated from the 
atomic density matrix, by summing the force contributed by each individual coherence. 
For the 2-(4)-2 level system, this is found using the following equation; 

W-J 2 2 y/3 22/ 22 /V<J 2 2 

-27I^^r-^-2^^r^) 

Calculations show that this system does not possess a cooling force in the steady state. 
The absence of force is a result of coherent population trapping, leading to formation 
of a velocity-independent dark state which inhibits the cooling cycle. However, time 
dependent calculations show that a transient force is present, as seen in figure 6, and 
decays to zero as the steady state is approached. Although the ground-state populations 
are found to continue evolving past 7t = 500, this transient force has decayed to zero by 
around 7t = 40, corresponding to the peak in gi/2 and G3/2 populations in figure 5(d). 
The initial force observed (at 7t = 1) is found to continue oscillating at larger velocities, 
with a regular frequency but low amplitude, before dying away. 




Figure 6. Time dependent dipole radiation force for the 2-(4)-2 multilevel atom 
(figure 2) as a function of normalized velocity. Saturation parameter G = 1 and 
detunings Si = —87, 82 = —27. 
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6. Optimal peirameters 

Analysis of individual coherences in the atomic density matrix shows that, for this 
system, the majority of the force is contributed by the ^Pi/2 — S1/2 transition. 
Calculations were also carried out to simulate cooling if all population is transferred 
to the ^P3/2 state, for example by a STIRAP process. These conditions led to a much 
reduced initial force. This indicates that to maximize the force experienced by the atoms 
it is preferable to have as large a proportion of the population as possible in the ^Pi/2 
state. This would need to be taken into consideration in any experimental proposal. 

Laser cooling schemes provide an extremely large parameter space, and as an 
example, we have investigated the laser detuning, 6. Figure 7 shows how the maximum 
initial force experienced by the atoms varies with laser detuning. The force is measured 
at resonant velocity, kv = Si = 62- It can clearly be seen that increasing the detuning 
above —27 does not result in a significant increase in force. For symmetry reasons the 
force decays to zero at kv — 6 — 0. 




Figure 7. Variation of the maximum radiation force with laser detuning, at resonance 
velocity, kv = 5; Saturation parameter G = 1, detuning 6 = di = 62 

It has also been observed that a small detuning difference between the two lasers 
is beneficial. The force is found to decay more slowly when the two detunings vary by 
a small amount, approximately when {61 — 62) = — O.57. It is possible that this small 
detuning difference minimizes interference effects between the two laser beams. 



7. Conclusion and Outlook 



An improved method for calculating the transverse cooling of atoms with complex 
groTind-statc structures has been presented, and the optimal detunings determined to 
maximize the force on an atom such as ^^Ga. The method is efficient and can be adapted 
to considerably larger multilevel systems. As an example of the kind of complex system 
that can be considered, the A-shaped cooling cycle outlined above is applied to a Group 
13 isotope with nuclear spin / = 3/2, for example ^^Ga. The ground state and first 
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excited state contain 32 magnetic sublevels, requiring 1024 density matrix equations for 
a full treatment. Using the method presented this calculation requires the formation 
of only three 32 x 32 matrices. One suggested cooling cycle is a 3-5-(l-3-5)-3 scheme, 
involving only one of the hyperfine excited states. This requires five different coloured 
lasers, as shown in figure 8(a). The evolution of ground-state populations at resonance 
velocity kv = 6 are shown in figure 8(b). For clarity only a selection of states are shown, 
and are labelled according to quantum numbers \J,F,Mp). 




Figure 8. (a) Energy level schematic for ^^Ga. Transitions forming a possible 3- 
5-(l-3-5)-3 five-colour cooling cycle are highlighted, (b) Example of the evolution of 
selected ^Pi/2, ^P3/2 ground-state populations in ^^Ga, at resonance velocity kv = S. 
Assumptions: all laser saturation parameters Gi — 1, all radiative widths are equal 
7i = 1, and all detunings Si = —2"/^. For clarity only a selection of states are shown, 
and are labelled according to quantum numbers \J, F, Mp) . Initial populations are 
distributed equally amongst magnetic sublevels in each fine structure ground state, 
with 80% in the state, and 20% in the ^P3/2. 



An interesting atomic system for future investigation is that found in the stable 
isotopes of Group 13 atom thallium, ^os^jj g^^^ ^°^T1, which both have 1=1/2 and are 
important candidates for sensitive measurements of parity non- conservation [24]. 
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In principle this method could also be used to model the laser cooling of molecular 
systems in which a suitable cooling scheme has been identified. The simplicity of the 
method means that electric or magnetic fields can be incorporated into the Hamiltonian 
matrix with little difficulty, as can time-dependent pulsed or shaped laser fields, as well 
as complex repumping schemes. This could be of great use in investigating methods for 
overcoming dark state resonances which hamper the laser cooling process. 
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